
static const float pi = 3.14159265f;



float3 ray_sphere_intersection(float3 c, float r, float3 o, float3 d)
{
	float3 res;
	
	float3 l = c - o;
	float s = dot(l, d);

	float l_2 = dot(l,l);
	float r_2 = r*r;


	if(s < 0.0f && l_2 > r_2)
		return float3(0.0f, 0.0f, 0.0f);

	float m_2 = l_2 - s*s;

	if (m_2 > r_2) 
		return float3(0.0f, 0.0f, 0.0f);

	float q = sqrt(r*r - m_2);

	float t;

	if (l_2 > r_2) 
		t = s - q;

	else
		t = s + q;

	return o + t*d;

}


float border;
float globe_radius;
bool bEnd;


float fog_factor(float3 ws_pos, float p_ws_y, float h, float3 vEye)
{
	float3 cam = vEye;
	float3 p = ws_pos;
	
	// fog should stop at sea level
	if (p.y < -0.5) return 0.0f;

	// if outside sphere...
	if (length(cam) > 5.0f)
	{
		cam = ray_sphere_intersection(float3(0.0f, 0.0f, 0.0f), 10.0f, cam, normalize(p - cam));
	}

	// modify max fog height
	h += cos( ws_pos.x - border * 5.5f/*multiplier works as wind speed*/) * sin(border) * 0.1f/*smooths out the cosine function*/;		
	h = clamp(h, 2.9, 3.1);

	float cam_y = min(cam.y, h);
	float p_y	 = min(p.y,  h);
	float deltaY = abs(cam_y - p_y);
	float deltaD = sqrt(pow(p.x - vEye.x, 2) + pow(p.z - vEye.z, 2));
	float dI = deltaY * deltaY * 0.5f; 
	float l = length(float3(-1.0f, 0.0f, 1.0f) - p);
	float xx = min(5.0f, ws_pos.x - border) + 5.0f;
	float att = 0.0f;
		
	// dense cloud
	att = 1.0f / ((xx * 0.5f * xx) + 1.0f);
	
	// dense/sparse cloud
	float k = max(att, 0.0f) + 
				(1.0f + 
				// sin function for variation
				max(sin(0.8f * pi * ws_pos.y), 0.0f) * 0.5f  + 
				// 2nd sin function, animated density (scaled by cos(border))
				max(sin(2.0f * pi * ws_pos.y), 0.0f) * 0.5f * max(cos(border), 0.3f)
				// sin function sum is scaled so that dense cloud effect is visible
				) * 0.4f;
	dI *= k;

	float density;

	if (deltaY - 0.001f < 0.0f)
		density = (cam_y < h) ? sqrt(1.0f + pow(deltaD, 2)) * k : 0.0f;
	else
	{
		density = sqrt(1.0f + pow(deltaD / deltaY, 2)) * dI;

		// fix for deltaY -> 0
		if (cam_y < h && p_y < h)
			density = max(density, 0.2f);
	}
	
	// rolling in
	
	// ( assume maximum globe radius (=10) )
	// ws_pos.x : [-10, 10] -> [0, 20]
	float r = ws_pos.x + 10.0f;
		
	float r_end = 0.0f + border;	
	float r_start = r_end - 8.0f;

	float falloff;


	//         1.0f ->       r_start         [varying]        r_end      <- 0.0f
	// -------------------------|-------------------------------|---------------------	
	if (r <= r_start)
		falloff = 1.0f;
	else if (r >= r_end)
		falloff = 0.0f;
	else
		falloff = (r_end - r) / (r_end - r_start);

	return density * falloff;
}


float4 computeFog(float4 in_color, float3 ws_pos, float ws_pos_y, float h, float3 vEye)
{
	float4 color;

	float fog_f = exp(
						-fog_factor(ws_pos, ws_pos_y, h, vEye)
					 );

	float4 fog_color = float4(0.9f, 0.9f, 0.9f, 1.0f);

	if (border > 20.0f)
	{
		float p = min((border - 20.0f) / (5.0f), 1.0f);
		if (p == 1.0f)
			return in_color;

		fog_f = fog_f * (1.0f - p) + 1.0f * p;
	}

	color = fog_f * in_color + (1.0f-fog_f) * fog_color;
	color.w = in_color.w;

	return color;
}

